Collective mode of homogeneous superfluid Fermi gases in the BEC-BCS crossover 

R. Combescot a , M.Yu. Kagan b and S. Stringari c 
(a) Laboratoire de Physique Statistique, Ecole Normale Superieure, 24 rue Lhomond, 75231 Paris Cedex 05, France 

and Institut Universitaire de France 
(b) P.L. Kapitza Institute for Physical Problems, Kosygin street 2, Moscow, Russia, 119334 
(c) Dipartimento di Fisica, Universita di Trento and BEC-INFM, 1-38050 Povo, Italy 

(Received February 6, 2008) 

PACS numbers : 05.30.Fk, 32.80.Pj, 47.37.+q, 67.40.Hf 

We perform a detailed study of the collective mode across the whole BEC-BCS crossover in 
fermionic gases at zero temperature, covering the whole range of energy beyond the linear regime. 
This is done on the basis of the dynamical BCS model. We recover first the results of the linear regime 
in a simple form. Then specific attention is payed to the non linear part of the dispersion relation 
and its interplay with the continuum of single fermionic excitations. In particular we consider in 
detail the merging of collective mode into the continuum of single fermionic excitations. This occurs 
not only on the BCS side of the crossover, but also slightly beyond unitarity on the BEC side. 
Another remarkable feature is the very linear behaviour of the dispersion relation in the vicinity of 
unitarity almost up to merging with the continuum. Finally, while on the BEC side the mode is 
quite analogous to the Bogoliubov mode, a difference appear at high wavevectors. On the basis of 
our results we determine the Landau critical velocity in the BEC-BCS crossover which is found to 
be largest close to unitarity. Our investigation has revealed interesting qualitative features which 
would deserve experimental exploration as well as further theoretical studies by more sophisticated 
means. 



I. INTRODUCTION 

The field of ultracold fermionic atoms has seen recently a remarkable achievement. With 6 Li as well as with 40 K, it 
has proved to be possible to display an experimental realization of a system going continuously from a Bose-Einstein 
superfluid to a BCS one [1]. This has been done by taking advantage of a Feshbach resonance, which allows to control 
by a magnetic field the scattering length, and thereby the effective interaction, for two fermionic atoms belonging to 
different hyperfine states (hereafter referred to as spin up and spin down states, as it is commonly done by analogy 
with electrons forming pairs in a superconductor). On one side of the resonance the scattering length a is positive and 
the spin up and spin down atoms can form diatomic molecules. At low enough temperature these bosonic molecules 
can undergo a Bose-Einstein condensation which has indeed been observed for 40 K [2] as well as for 6 Li [3-5]. On 
the other side of the resonance where the scattering length is negative, evidence for BCS superfluidity has also been 
observed, for example through pair-breaking features in collective modes [6], observation of the pairing gap [7] and 
more directly and recently observation of vortices [8]. When the magnetic field is tuned across the resonance, where the 
scattering length diverges, the superfluid undergoes a crossover from Bose-Einstein condensates to a BCS superfluid. 
Around resonance the superfluid is, so to speak, of a new kind since it is neither a full BCS superfluid nor a full 
molecular condensate. It is obviously quite interesting to study in detail this new superfluid, although it is clear that 
it will share the properties common to a BEC and a BCS superfluid. 

In this paper we will be interested specifically by the way in which the Bogoliubov collective mode, which is well 
known for Bose condensates since it is the only elementary excitation allowed in this system, goes over into the 
so-called Anderson-Bogoliubov mode which is known to exist in the BCS case. In the limit of very low frequency 
this is easy to understand since in both limits this mode is just identical physically to sound waves, propagating in 
the neutral superfluid, which are coupled to fluctuations of the phase of the order parameter. Naturally the same 
interpretation will also hold in the crossover. However we will be interested in this paper with all the range of possible 
frequencies, which will clearly display a larger number of interesting features. In particular it is worthwhile to enquire 
how the composite nature of the boson corresponding to the mode will show up in the dispersion relation. 

Our study of the collective mode will rely on the use of dynamical BCS theory. Naturally this theory is not exact 
in the crossover, while it is known for quite a long time to give a correct description not only in the BCS limit, but 
also in the deep BEC regime [9]. In between it should be regarded as an interpolation model, which has no specific 
reason to be quantitatively reliable. On the other hand BCS theory is a very coherent framework which is known to 
display all the qualitative physics on both sides of the crossover and it is therefore reasonable to trust it qualitatively 
in the crossover, and to make use of it at least as a first instructive step. Actually it has already been used in the 
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literature for some partial studies of the collective mode. For example Minguzzi, Ferrari and Castin [10], calculating 
the dynamic structure factor, have stressed that the observation of this mode in a trapped ultracold gas would be 
a clear signature of superfluidity. BCS theory has also been used by Biichler, Zoller and Zwerger [11] to consider 
essentially the low frequency regime. Finally Pieri, Pisani and Strinati [12] had also to consider the collective mode, 
in the course of their calculation of the spectral function which was their major focus and which was mostly aimed 
at the physics of high T c superconductors, but they did not proceed to a detailed study as we will do in the present 
paper. 

The paper is organized as follows. In the next section, as an introduction, we recall for the static properties the 
standard BCS model widely used to describe the BEC-BCS crossover. We take advantage of it to introduce some 
new results concerning the sound velocity and also some simple relations which appear within this model around the 
yti = point. Then we give the principles behind the proper handling of the BCS dynamics, the technical details being 
deferred to the appendix. We present the result for the collective mode dispersion relation and we give our results for 
the detailed study of this mode, on the BCS side, on the BEC side and around unitarity. Once again, for clarity and 
simplicity, the technical analytical details related to this study are given only in the appendix. 



II. STATIC BCS FRAMEWORK 



Let us first briefly recall the standard BCS framework and how it is applied to the BEC-BCS crossover, in the way 
initiated by Leggett [13] and by Nozieres and Schmitt-Rink [14], and devclopped in particular by Sa de Melo, Randeria 
and Engelbrecht [15,16]. This BCS model has been for example used recently [17] by Viverit et al, together with a 
semiclassical approximation, to calculate the momentum distribution in trapped gases across the BEC-BCS crossover. 
This is often called the mean-field BCS model. However there is no diagonal mean-field included in this framework 
[18], the only mean-field character being just the fact that BCS theory itself is formally a mean-field theory. For this 
reason we drop the "mean-field" qualificative since it is redundant at best, and misleading at worst. 

In the BCS framework one eliminates the usual high energy cut-off oj c and the attractive interaction V, in the 
standard gap equation: 

in favor of the scattering length a, which leads to: 

— =Y(- —\ (2) 
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where e k = k 2 /2m is the atomic kinetic energy and E k = \f£%. + A 2 the single particle excitation energy. Here we 
have set = e k — /i, with fi being the chemical potential and A the gap parameter, and we take K = 1 throughout 
the paper. This gap equation is supplemented by the equation giving the atom density n for a single spin state: 



This equation, which does not play any role in the standard BCS weak coupling theory, is here an essential ingredient. 
Together with Eq.(2), it gives the evolution of the gap and the chemical potential, as a function of l/k F a, across the 
BEC-BCS crossover. After some manipulations we rewrite Eq.(3) as: 



A 2 f°° k 4 

k 3 F = — Ja Ja= dk — ^ (4) 
2m J Q E k 3 



Correspondingly we can rewrite Eq.(2) as: 
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In the BEC limiting case where [i — > —00 and A/fj, — ► one has J 2 = (7r/16).(2ra/|/z|) 3 / 2 and J 4 = 
(37r/16).(2m) 5 / 2 /| M | 1 / 2 . One finds naturally from Eqs.(2-3) H = l/(2ma 2 ), and A/\n\ = (32nna 3 ) 1 /' 
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A. Vicinity of n — 



It is of interest to consider more closely the vicinity of /j = 0, which is obtained for l/(kpa) — 0.553, since it is 
around this region that the gas will switch from a typical BCS behaviour to a typical Bose one. We consider the 
results of first order expansion in the vinicity of /x = 0. Details are given in the Appendix. It turns out that all the 
specific numerical integrals coming in this calculation can be expressed in terms of Euler's gamma function T(l/4). 
This leads to very simple relations between various physical quantities in terms of the universal parameter l/(kpa)o 
or equivalently (A/Ep)o- One finds indeed that the gap to Fermi energy ratio satisfies merely: 
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where the index indicates that the relevant physical quantities are taken for /x = 0. The behaviour of A/Ef and 
1/fcpa in the vicinity of \i = 0, can be simply expressed in terms of (A/Ep)o or, equivalently (l/fcF«)o- One finds: 
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It would be interesting to compare these results with the ones of Monte-Carlo calculations and with experiments to 
check how good is the quantitative description of BCS theory in this respect. 



B. Sound velocity 

Finally let us consider in this model the sound velocity c s , given by c 2 — (n/m)(dfi/dn). Varying Eq.(2) gives 
d/i/dA = J 2 A/J^, where we have introduced: 

J, = A-J 4 -»J 2 = J o dk k -^-rnJ Q dk i- (9) 

while varying Eq.(3) leads to 4ir 2 (dn/d(j,)A = J2A 2 and 4tt 2 (dn / 'd A) ^ = J^A. Taken together with Eq.(4) this leads 
to: 

2 _ 1 J 2 J4A 2 
Cs ~ 3m 2 J 2 A 2 + J 2 [ ' 

In the BCS limit, where fi ~ Ep = k 2 F /2m > and A/it — > + , one finds from Eq.(4) J4 = k 2 F Ji = 2rak\j A 2 while 
~ — (2m 2 /kp) ln(A/Ep) is negligible due to particle-hole symmetry, leading to the well known result c 2 = kp/3m 2 . 
On the other hand the dependence of n on A, i.e. J^, can not be omitted in the rest of the BEC-BCS transition. 
For example at unitarity, Eq.(5) implies J 4 = 2mJ 2 (A 2 + tt 2 )/it, while accordingly Eq.(9) gives = J2A 2 / 'fi. 
When this is inserted in Eq.(10) this merely leads to c 2 = 2/x/(3to). But this is what is expected at unitarity 
from universality arguments proving that, for dimensional reasons, \i is proportional to kp, i.e. n 2 / 3 , which give 
c 2 = (n/m)(d^/dn) = (2/3)/i/m. Hence, at unitarity, BCS theory leads naturally to a result in agreement with 
universality, but only provided that the important dependence of n and /x on A is retained. 

This is also the case in the BEC limit, where A//i — > 0_. Indeed we have seen that J 4 = Qm\fi\ J2 and = 4|/x| J2, 
leading to c 2 = A 2 /(8m|^|) = 2ima/m 2 . This coincides with the standard sound velocity c 2 = gunu /tum of a 
molecular BEC gas with % = n, vtlm — 2m and gM — faraM l m M where om = 2 a is the value for the molecular 
scattering length predicted by the BCS model [15], as opposed to the exact result [19,20] «m = 0.6 a. 
Finally in the case \x = we end up with: 



k 2 F (tt/2)(1 / 'kpa) 2 „ nl , fc 2 
s 2,m 2 l + ^ 2 /%){l/k F a)l~ m 2 

which can also be obtained from the relationship c 2 = (n/m)(d^/dn) with /j,(n) given by the expansion Eq.(8). 
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III. BCS DYNAMICS 



In order to obtain the collective mode, we have now to deal with the dynamics of the BCS theory. This can be done 
following various approaches, which are actually all basically equivalent presentations of the standard self-consistent 
BCS theory. To be specific we will consider here the kinetic equation approach. More details on this method will 
be given in the Appendix, where we will also sketch the derivation following the Green's function, or diagrammatic, 
approach to the BCS theory. 

The important physical point about this dynamics is that it has to be self-consistent. This means that not only 
the single particle excitations are set in motion, but the order parameter has also space and time dependences. This 
already appears clearly from the results given in the preceding section for sound velocity where, in order to obtain 
the proper result, one has to consider (dn/dA)^, in addition to (dn/dfi)A- This shows that one has to consider 
fluctuations of the modulus |A| of the order parameter. Similarly one has also to consider fluctuations of the phase 
of the order parameter, since its gradient is directly linked to the superfluid velocity. Hence the proper physical 
description is the dynamical evolution of the excitations in the presence of a dynamical order parameter, which is 
itself self-consistently determined by the excitations. If we did not take these fluctuations of the order parameter 
into account, we would have [21] a dynamics which would violate gauge invariance and particle conservation, and at 
the same time would not display the gapless collective mode. It is well known [21] that the existence of the gapless 
collective mode is intimately related to gauge invariance. 

These features are quite clear in the kinetic equation formulation of the BCS theory, which is a generalization to the 
superfluid case of the standard Landau's kinetic equation for a Fermi liquid. In addition to the diagonal self-consistent 
fields present in Landau theory, superfluidity introduces off-diagonal self-consistent fields and leads to deal with the 
following 2x2 matrix for the single particle density: 



where A, ji = 1 or 2, and by definition c k ,i = Cfe-f and c k , 2 = c _fei are Nambu doublets. At the level of small 
fluctuations with wavevector q, the fluctuating part 6n k of this density matrix satisfies the kinetic equation: 



which is formally identical to the standard kinetic equation, except that one deals now with 2x2 matrices. The 
equilibrium energy matrix e k , together with its fluctuating part 5e k , and the corresponding equilibrium density matrix 
n£ are given explicitly in the appendix. 

As usual we find the dispersion relation for the collective mode by looking at the momentum and frequency (q, u) 
which makes the response of the system infinite, as for any eigenmode. With respect to the specific response function 
to be considered, one may think of calculating the density-density response function since in the low frequency limit 
the collective mode corresponds physically to sound waves [citer Randeria et travaux antrieurs]. Or else one may 
consider the anomalous response function since, on the BEC side for example, the collective mode coincides with the 
Bogoliubov excitations which are just the elementary excitations of the system out of the condensate. Actually both 
response functions lead naturally to the same answer for the collective mode, because density fluctuations are coupled 
to phase fluctuations of the order parameter since the particle current is proportional to the gradient of the phase. 
Naturally the response function is by itself quite interesting, not the least since its imaginary part is directly related 
to the dynamical structure factor, which has been studied very recently in the BEC-BCS crossover. We give also its 
explicit expression in the appendix for completeness. 

In the following section we give the explicit expressions allowing to calculate numerically and study analytically 
the collective mode dispersion relation. Some details on its derivation are given in the appendix. We then go on and 
study it in various situations along the BEC-BCS crossover. 



From the kinetic equation or the diagrammatic method (see Appendix) we end up with the following dispersion 
relation, between the frequency ui of the collective mode and its wavevector q: 



(12) 
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IV. COLLECTIVE MODE DISPERSION RELATION 



T T 2 r2 



(14) 



where the various integrals / arc defined by : 
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and we have set E± = + A 2 , E = ^C 2 + A 2 with £± = (k 2 ± fcgu + <7 2 /4)/(2to) - fi and £ = k 2 /(2m) - fj,. In 
the above integrals we have introduced the variable u = cos 9 where 9 is the angle between k and q. 
Now it is quite interesting to note that one can also write: 



J n (w,g) = i f°°k 2 dk f du { ^±- 
1 Jo Jo 
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The identity between Eq.(16) and Eq.(18) does not come by chance, but it is just a consequence of mass conservation 
which is naturally satisfied by our starting equations. We could have used this mass conservation in the course of 
our derivation to find directly Eq.(18), but it is here more efficient to just check explicitly the identity. Indeed this 
identity is equivalent to state that: 



E 
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(E++E-) 2 -uj 2 



— - V-f— — - — 

E~2-^2^E: + E-' E 







(19) 



where the last equality is obtained by making appropriately the changes k — > k ± q/4 in the summation variable. 

Result (18) for Iu(u>, q) makes it very easy to check that the dispersion relation (14) exhibits the proper phononic 
behaviour ui = c s q as uo and q tend to 0. Indeed it this case we have /12(C), 0) = (1/2) J 5 , ^22(0, 0) = — A 2 J 2 and 
Iu(oj, q) ~ oj 2 J2/4 — q 2 J4/ (12m 2 ) which leads immediately to the expression Eq.(10) for the sound velocity. 



V. RESULTS 



A. THRESHOLD FOR SINGLE FERMION EXCITATION SPECTRUM 



All the above integrals display singularities when lu = E + + E-. This equality corresponds physically to the 
possibility of breaking a Cooper pair into two fermionic excitations. The spectrum for this kind of excitation is 
continuous. The corresponding domain is bounded from below by the line u>th = m in (E+ + E_) where the minimum 
is to be taken over all the possible values of k. It is easy to see that the minimum is obtained for k.q = which 
gives E+ = E-. In the case yu > and for q < 2^/2m^, the minimum of E + + E- is reached for k 2 /2m = 
fi — q 2 /8m, so that co t h = 2A. For fi > and q > 2yJ2m[i, the minimum is reached for k = 0, which leads to 
w t h = (q 2 /8m — [i) 2 + A 2 . On the other hand, for the case fi < 0, the minimum is always reached for k = 0, 
leading to u t h = 2\/ (q 2 /8m + \fi\) 2 + A 2 . To summarize we have: 

wth = 2A for (i>0 and ?<2 y^2mp (20) 

LOth = 2\/ (q 2 /8m — n) 2 + A 2 otherwise 

A few examples of this threshold are given in Fig.l. This result differs from the predictions of the non- interacting 
Fermi gas result, which has not only a lower boundary u>th- = for q < 2kp, and coth- = —qkp/m + q 2 /2m for 
q > 2kF, but also an upper boundary uJth+ — qkp/m + q 2 /2m. It is worth noticing that this behaviour is not recovered 
in the A — > limit reflecting the fact that, with respect to the excitation spectrum domain, the normal state limit is 
singular. However one recovers of course in this limit the correct behaviour of the dynamic structure factor, because 
outside the standard normal state domain the structure factor will go to zero, due to the effect of coherence factors. 

In the BEC limit, one has A/fi — > 0_ which gives to t h — 2|/x| + q 2 /Am. The physical meaning is clear. This is 
just the energy necessary to break a molecule with binding energy = 2|/j| = 1/ma 2 into two fermions with total 
momentum q, leading to an additional kinetic energy (l/2)g 2 /(2m), just as for a molecule of mass 2m. 
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FIG. 1. Threshold uj th as a function of wavevector q for, from bottom to top, 1/fcjro = -1. , -0.5 , 0. (unitarity) , 0.553 
((i = 0) , 1. and 2. 



B. WEAK COUPLING BCS REGIME 



Let us consider first the weak coupling regime on the BCS side, where a — > 0_ and A — > 0. The numerical result 
for 1/kpa = — 1. and —0.5 are displayed on Fig. 2. Although in these cases a is not so small, they still display a 
behaviour which is typical of this limit and which would just be more strongly marked at smaller a. For increasing 
q the dispersion relation has first a fairly extended linear part, the slope being naturally the sound velocity. Then it 
bends to become nearly horizontal. This last behaviour is quite reminiscent of an anticrossing in a two level system. 
Here the collective mode plays the role of one of the levels. The other one is not really a level since it is rather the 
pair-breaking continuum. However the well-known divergence of the BCS density of states at the uo = 2A threshold 
makes the threshold level quite similar to a single level at energy u = 2A. It is quite interesting to note that this 
behaviour of the dispersion relation of the collective mode is a unique feature of the ultracold Fermi gas system. Indeed 
in superconductors the collective mode frequency is pushed up to the plasma frequency because the paired electrons 
are charged particles. On the other hand, in superfluid 3 He which is a neutral BCS superfluid, the strong hard core 
repulsion between atoms produce a much faster sound velocity and the mode actually merges in the pair-breaking 
continuum. Hence the anticrossing feature is found only in Fermi gases. 

Looking at Fig. 2 one would wrongly conclude that the collective mode touches the pair-breaking continuum and 
merges with it for q < 2^/2mJl. This is not the case, as the above anticrossing argument makes it clear and can be 
analytically shown (see the Appendix). As we will see in the next section V C, the collective mode does indeed merge 
in the continuum in the weak coupling regime, but this always occurs for q > 2^/2rnjl, the merging point approaching 
the q = 2y/2m(i point as the coupling becomes weaker and weaker. 



6 




FIG. 2. Dispersion relation uj/Ef of the collective mode as a function of q/kF for 1/kpa = — 1. (lower thick line) and —0.5 
(upper thick line) . The location (see section V A) of the threshold for pair- breaking is given in each case by the thiner line. 



C. CONTACT OF THE COLLECTIVE MODE WITH THE CONTINUUM 



In order to clarify the position of the collective mode dispersion relation with respect to the threshold of the 
continuous excitation spectrum, it is useful to find precisely the wavevector q m (and corresponding frequency w m ) at 
which they meet. Since the meeting never occurs for \i < one can set to = 2yJ (q 2 /8m — it) 2 + A 2 with it > in 
Eq.(14) and solve the corresponding equation numerically The result is displayed on Fig. 3. 

First of all, it is clear from these results that, as soon as a < 0, the mode meets the continuum slightly beyond the 
point u) = 2A, and q = 2^/2mit, where the lower boundary of the continuum begins to depart from the simple value 
cj t h = 2A. The corresponding wavevector is indicated as a dashed line on Fig. 3. When one goes towards large and 
negative values of a, the meeting point (q m , u m ) gets extremely close to (2^2m/j,, 2A), while never reaching it exactly. 
The meeting point begins to depart appreciably from this location in the q — lo plane only when one gets close to 
unitarity. However at this stage a rather surprising feature appears. On the other side of the resonance, for a > 0, 
a second branch of the dispersion relation for the collective mode takes place at very large values of q, extending up 
to q — > co. The appearance of this large q branch at unitarity can be confirmed analytically. This branch meets the 
continuum at some lower value q m . This means that the above equation has two solutions instead of one, as it was 
the case for a < 0. This second solution corresponds to the upper branch of the curve seen in Fig. 3. When one 
goes away from unitarity on the a > side, the two meeting points move rapidly towards each other. They are found 
numerically to merge into a double solution q/kp — 2.57 for 1/kpa — 0.161. 

Beyond this value of 1/fcpa, the above equation has no longer solutions, which means that the collective mode 
never meets the continuum. Hence, qualitatively, one has already the situation found far on the BEC side. It is 
remarkable that this evolution occurs very rapidly, as a function of the scattering length, since this happens basically 
between 1/kpa = and 1/kpa — 0.16. It is also quite noticeable that this change occurs while the chemical potential 
is still positive, a regime where one would not be tempted to speak of molecules since the molecular formation would 
seemingly imply a negative chemical potential. Nevertheless in this regime we find a fully characterized collective 
mode, similar to the elementary excitations of the molecular condensate. 
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-0.2 0.2 

1/k F a 

FIG. 3. Wavevector q m , in units of kp, at which the collective mode dispersion relation meets the excitation continuum, as 
a function of 1/k.Fa (note that the q m scale is logarithmic). The dashed line corresponds to q = 2y / 2m7i. It goes to qjkp = 2 
when a —t — oo 

It is quite tempting to interpret the remarkable fact that the collective branch meets the continuum twice for large 
and positive values of the scattering length a in the following way. It has been shown [26,27] that, in the normal state, 
the presence of the Fermi sea of atoms makes it more difficult for molecules to form, because some states are blocked 
due to Pauli principle by the presence of other atoms. This is mostly felt for molecules with low total momentum, 
which do not form at unitarity, but rather only at some positive value of a depending on temperature. On the other 
hand molecules with very high momentum are insensitive to the presence of the Fermi sea and form right at unitarity. 
Now, for a standard Bose-Einstein condensate, the physical nature of the collective mode at very high q is just the 
excitation of a single boson of momentum q out of the condensate. Similarly for the Bose condensate of fermionic 
molecules occuring in our full BEC regime, the collective mode at very high q is the excitation out of the condensate 
of a molecule of momentum q. The natural extension of the above result, regarding molecular formation in the normal 
state, is to say that, just for the same physical reasons, these very high q excited molecules form also unhindered in 
the superfluid state at T = 0. Hence for very high q the collective mode exists as long as these molecules exist, that 
is up to unitarity when we increase a. In other words for positive values of a and large q one always expect to find a 
discretized branch due to the fact that in this case the two-body problem always admits the existence of a molecular 
state and many body effects are negligible. Naturally the difference between this excitation and molecular breaking, 
which corresponds to the lower boundary of the continuum, decreases when a increases and at unitarity, where the 
molecular binding energy becomes zero, the collective mode merges into the continuum. On the other hand when we 
decrease q, we expect that the formation of excited molecules will become more difficult due to Pauli exclusion by the 
condensate, because the Pauli principle becomes more and more effective, reducing effectively in this way the binding 
energy. This explains why at some stage the collective mode merges into the continuum when q is decreased. Let us 
stress that, if this physical interpretation is correct, this implies immediately that the exact theory of the BEC-BCS 
crossover should also display this behaviour because our arguments do not make use specifically of the BCS theory. 
This means that we exclude a scenario where the full curve in Fig. 3 goes to infinity as one approaches the resonance 
from the negative a side, which would result in this curve being entirely located on the a < side. Clearly it would 
be worthwhile to investigate experimentally this question. 

D. BEC LIMIT AND HIGH FREQUENCY LIMIT 

We expect naturally to find the Bogoliubov mode when we go to the BEC limit, which corresponds to fj, < and 
a — > 0+. Physically we have in this case a dilute gas of very tight molecules, which should behave as elementary 



bosons. In this limit the well-known dispersion relation of the Bogoliubov mode is given by: 



2 2,/' <? 2 
c s q + 



V2M„ 



(21) 



where, as we have seen, the sound velocity c s is given by c 2 — g m n m /M m — 2ima/m 2 and M m = 2m is the molecular 
mass. Since we are in the a — > 0+ limit, the first term will be negligible compared to the second one. More precisely 
this means that we will consider wavevectors q 3> (na) 1 / 2 , i.e. much larger than the inverse of the healing length 
of the molecular Bose gas. In this case we will have merely w ~ q 2 /Am and we are in the regime where the mode 
corresponds physically to kick a molecule out of the condensate. 

In this BEC limit the gas parameter na 3 becomes very small and A/|/i| ~ (na 3 ) 1 / 2 — > 0. On the other hand the 
above condition q 3> (na) 1 / 2 is equivalent to q 2 /2m 3> A. Hence, in this limit A becomes negligible compared to all 
the other energies in the problem, which means that we have to take it formally going to zero. In this situation, we 
have merely to take E± = £± and the integrals in Eq.(14) can all be easily performed analytically. One finds: 



Ji2 = {2mf 2 ^- \ V%| + q 2 /2m + 2lu - ^ A\^i\ + q 2 /2m - 2oj 



(22) 



hi = I22 = (2m) 



3/2 



1/2 _ 



+ q 2 /2m + 2u> - y/4\n\ + q 2 /2m - 2u> 



(23) 



Then Eq.(14) becomes hi +u)I\2 — (since hi < 0) and one finds immediately that the solution is indeed w = q 2 /Am. 

We note also that, if we investigate the regime where q and to go to infinity whatever the value of na 3 (provided 
fj, < 0), we are lead to perform exactly the same calculation. Indeed in this limit, u> and q 2 /2m are large compared 
to A, so we can take again A — > 0. Therefore the result for the dispersion relation is also lo = q 2 /Am. This is 
again in agreement with the physical picture holding for the mode in this limit, namely kicking a molecule out of the 
condensate. 



E. THE COLLECTIVE MODE AT UNITARITY AND ON THE BEC SIDE OF THE CROSSOVER 

Now we will proceed to give results for the collective mode dispersion relation on the BEC side of the crossover. 
We first exhibit the result right at unitarity (see Fig. 4) . It is somewhat analogous to what has already been found on 
the weak coupling side in section VB, with the dispersion relation merging into the continuum for q/kp = 1.76. 

A noticeable feature of Fig. 4 is that the dispersion relation is surprisingly almost linear up to the merging with the 
continuum. This result predicted by the BCS model can be relevant for understanding the experimental data on the 
collective oscillations in trapped Fermi gases. In fact in actual experiments carried out in elongated traps at unitarity 
the superfluid gap A is never much larger than the radial oscillator frequency lu±_ which, within hydrodymamics, fixes 
the collective frequency of the radial breathing mode according to the law [23] uj ra d = y/l0/3u>±. This means that in 
these experiments one is exploring a region of relatively high wavevectors. The quasi-linearity up to high wavevectors 
allows one to understand why hydrodynamics (valid in principle at low wavevectors) is in such a good agreement with 
experiments [24,25]. A careful study of the dependence of the collective frequency on the ratio A/ui ra d might actually 
provide new insight on the behaviour of the excitation spectrum and in particular check if the extended linearity of 
the dispersion law is a peculiarity of the BCS model or is a more general feature. 

On Fig. 5 we display the results when fi = 0, corresponding to 1/kpa = 0.553 as we have seen in section II. Since 
we are beyond the value 1/kpa = 0.161, the dispersion relation is now fully below the continuum as discussed in 
section VC. However for most of the range it is quite close to lower boundary of this continuum, and the situation 
is not drastically different from the one obtained at unitarity. This is even more so if we take into account that the 
spectral weight associated with the collective mode gets smaller when it goes near the continuum. On the other hand 
the shape of the dispersion relation is already quite similar to the familiar one given by the Bogoliubov result, already 
discussed in section VD. We have made this clear by plotting the result of the Bogoliubov dispersion relation in the 
same figure, with the sound velocity c s taken naturally for /x = 0. Although we are certainly not in the deep BEC 
regime where this formula is expected to be valid, we see that it is quite similar to the actual result, although it goes 
unappropriately in the continuum. 
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FIG. 4. Dispersion relation lu/Ef of the collective mode as a function of q/kF at unitarity (thick line). The location (see 
section V A) of the threshold for pair-breaking is given by the thiner line. The collective mode merges into the continuum for 
q/kF = 1.76, as it can be seen from Fig. 3 



The situation becomes quite different when we move slightly deeper in the BEC regime as it can be seen in Fig. 6, 
where we have plotted the results for \jkpa = 1.. We see that the collective mode has moved much further from the 
continuum. At the same time the actual result is quite closely approximated by the Bogoliubov formula. Hence it 
is not surprising that, when we go to l/fopa — 2., the Bogoliubov result becomes almost undistinguishable from our 
numerical result, as it can be seen in Fig. 8, which shows that we are already effectively quite deep in the BEC regime. 
This is coherent with the fact that the absolute value \fi\ of the chemical potential is already almost 4 times the Fermi 
energy Ep, as it can be seen from Ref. [15,16]. Noting that for the dilute limit, this chemical potential is half the 
free molecular binding energy \jj,\ = l/2ma 2 , leading to \[i\/Ep — l/{kpa) 2 — 4, one sees that one is indeed already 
quite near this situation. This excellent agreement between the BCS result and the Bogoliubov formula makes it 
worthwhile to investigate the collective mode dispersion relation in the dilute limit to next order in powers of A/|/z|, 
that is to first order in (A/|/i|) 2 . This means going one step further than what we have done in the preceding section 
VD. This turns out to be not so complicated analytically. Details of the derivation of this expansion are given in the 
appendix. One obtains: 

^ M ^IfzM! + (jq 2 (24) 



4 16+(<7a) 2 \ v 4m / 

up to first order in (A/|yu|) 2 <~ na 3 , where we have used = l/(2ma 2 ) to lowest order. This has to be compared 
with the result from the Bogoliubov formula Eq.(21). Taking the sound velocity c 2 = A 2 /(8to|/i|) in the BEC limit 
from section II, we see that, in the limit qa <C 1, this is in perfect agreement with Eq.(57) and in particular we recover 
the well known result 

q 2 Aima , . 

ijJ = -f h 25 

4m m 

which expresses the positive shift of the dispersion with respect to the free value q 2 /4m in terms of the scattering 
length. On the other hand, in the limit qa 1, the dispersion relation becomes instead: 

Q 2 47ma 

w = (26) 

4m m 

showing that in this case the shift has a negative sign. These results are in qualitative agreement with our numerical 
results seen in Fig. 6 and 8. Indeed these figures show that the agreement with the Boboliubov result (which is 
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indeed above the numerical result) is not so good for large qa, whereas it is remarkably good for lower qa (note that 
in the linear range qa <C A/|//| we have automatically agreement since the sound velocity is the same). This result is 
physically reasonable. Indeed, for q <C 1/a, one is testing the system at lengths which are larger than the molecular 
size a. Hence it is natural that one gets the behaviour of pure bosons and an agreement with Bogoliubov theory. On 
the other hand for g> 1/a one is testing lengths which are smaller than the molecular size, and the fact that we deal 
with composite bosons appear at this scale. It is therefore natural to find a disagreement with Bogoliubov theory 
which does not include this physics. 




FIG. 5. Dispersion relation lu/Ef of the collective mode (thick line) as a function of q/k,F for /i = 0, corresponding to 
1/kFd = 0.553. The location (see section VA) of the threshold for pair-breaking is given by the thiner line. The dashed line 
indicates the result obtained from the Bogoliubov formula for the collective mode for this value of /i. 



F. Critical velocity 

We have now basically all the informations required to calculate the critical velocity of the superfluid v c in the 
BEC-BCS crossover. Indeed, according to Landau's criterion, it is given by v c — Mm(u)(q)/q) where uj(q) is the 
energy of an elementary excitation. For this BCS superfluid we have two kind of excitations. First we have bosonic 
excitations corresponding to the collective mode we are studying. If the dispersion relation has an upward curvature, 
the minimum of u{q)/q) is obtained for q — > 0, which gives for v c the sound velocity v s calculated in section II B if we 
retain onl y this kin d of excitations. However we have also to take into account fermionic single particle excitations 
uj(k) = y^fe + A 2 , and we have to find again the minimum of co(k)/k. However we have just seen that the lower 
boundary for the pair-breaking continuum corresponds to the creation of two single particle excitations with same 
energies and wavevectors, since we have found E + = E_ in this case. Hence it is easily seen that, since the wavevector 
corresponding to the minimum will clearly be greater than y/2m\fi\, Min(o>(fc)/fc) for single particle excitations is 
identical to Min(oJt/ l (g)/g) for the pair-breaking continuum found in Fig.l (the correspondence is just obtained by 
multiplying energies and wavevectors by a factor 2). This minimum is easily found analytically to be given by 
([\/ A 2 + (i 2 — yu]/m) 1//2 , obtained for the single particle wavevector [4m 2 (A 2 + /z 2 )] 1 / 4 . Hence, taking into account 
both kind of excitations, the critical velocity is given by: 

v c = Min (v s , ([VA 2 +m 2 - f4/™) 1/2 ) (27) 

which is plotted in Fig. 7. The result displays a kink at its maximum, occuring when one switches from bosonic to 
fermionic excitations. It occurs very near unitarity, on the BEC side. It is worth noticing that, in terms of critical 
velocity, the strength of the superfluid is at its highest around unitarity and not on the BEC side, as one might 
perhaps think at first. 
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In the above determination of the critical velocity, we have considered that the collective mode dispersion relation 
is bending upward. However, while this is clearly so when the chemical potential is negative (see Fig. 8, 6 and 5), 
this is no longer true at unitarity (see Fig. 4) and on the BCS side (see Fig. 2). Nevertheless these last situations 
are irrelevant for our purpose since in these cases, the critical velocity is due to single particle excitations. But one 
may still wonder if, near unitarity, one does meet the case where some part at least of the collective mode dispersion 
relation would display downward bending, which could result in Mm(oj(q) / q) obtained not for q — > 0, but for some 
finite value of q, as it occurs for example in superfluid 4 He for the roton minimum. To clear up this point, we have 
investigated in detail the collective mode dispersion relation in the vicinity of the maximum of v c . Indeed the switch 
between the behaviour found in Fig. 4 and in Fig. 5 occurs by an upward bending appearing in the low q region, while 
the higher q part has still a downward bending. However for this higher part we find always that ui{q)/q > c s so it 
does not play any role in the determination of the critical velocity, and accordingly there is no modification to bring 
to Fig. 7, and in particular the kink is not removed. It must also stressed that, for all this region near unitarity, these 
bondings (upward or downward) are extremely small, and the most striking feature is that the dispersion relation is 
remarkably essentially a straight line between q = and its lowest merging with the continuum. This is particularly 
so in the immediate vicinity of the maximum of v c where the difference with a straight line is quite minute and could 
not be seen on a figure (which accordingly spare us from displaying it). Nevertheless the numerical calculation shows 
that this is never exactly a straight line. It would clearly be very interesting to know if this peculiar and remarkable 
behaviour is just a coincidental specificity of the BCS model, or if it has a more general validity. 




FIG. 6. Dispersion relation lu/Ef of the collective mode as a function of q/kF for 1/ltFa = 1. (thick line), for which 
A/|/i| = 1.66. The threshold for pair-breaking is given by the thiner line. The dashed line indicates the result obtained from 
the Bogoliubov formula. 



An important physical quantity directly related to v c is the healing length, which we define merely as: 

? = — (28) 
mv c 

It has been studied in details by Pistolesi and Strinati [28]. We have plotted it in the insert of Fig. 7. It coincides, 
apart from a trivial numerical factor, with the usual definition of healing length in the molecular BEC limit and with 
the size of Cooper pairs, i.e. the coherence length, in the opposite BCS limit. The healing length, which takes its 
smallest value very near unitarity, fixes the core size of quantized vortices. 
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FIG. 7. Critical velocity v c (full line) and sound velocity c s (dashed line) as a function of 1//cfO- The insert gives = Vf/v c 



VI. CONCLUSION 



In this paper we have made a detailed study of the collective mode across the whole BEC-BCS crossover in fermionic 
gases at zero temperature. This has been done on the basis of the dynamical BCS model, which we have formulated 
equivalently in terms of diagrammatic theory or with collisionless kinetic equation. Naturally we have first recovered 
the results of the linear regime, where the mode is physically identical to sound propagation. We have cast the result 
for sound velocity in a simple form and payed particular attention to the vicinity of the point where the chemical 
potential is zero. Then we have turned to the non linear part of the dispersion relation. Our investigation has revealed 
interesting qualitative behaviours which would deserve further exploration with better theoretical tools, and would 
also be worthwhile to investigate experimentally. Quite generally there is an interplay between the collective mode 
and the continuum of single fermionic excitations. Roughly speaking the behaviour is analogous to the anticrossing 
phenomenon well-known for two level systems. This behaviour is most clearly displayed on the BCS side, which 
is quite interesting in itself since it provides the first experimental example of a pure Anderson-Bogoliubov mode. 
Nevertheless this anticrossing behaviour does not prevent the collective mode to merge at some stage into the single 
fermionic excitations continuum. This occurs mostly on the BCS side of the crossover, but quite interestingly it 
exists also slightly beyond unitarity on the BEC side, where we have made a careful study. This feature might be 
related to the physics of molecular formation in the crossover. Another remarkable feature, which is coherent with 
experiments on collective excitations in trapped gases, is the very linear behaviour of the dispersion relation in the 
vicinity of unitarity almost up to merging with the continuum. Finally on the BEC side the mode is naturally quite 
analogous to the Bogoliubov mode and no merging occurs anymore. However, even deep in this BEC regime, there 
is at high wavevectors a difference between this mode and the result from the BCS model, which can most likely be 
ascribed to the composite nature of the molecular bosons when they are tested at length scales comparable to their 
size. Finally our results on the collective mode and on the single particle excitations allows us to obtain the Landau 
critical velocity, which display a sharp maximum very near unitarity. 
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FIG. 8. Dispersion relation u)/Ef of the collective mode as a function of q/kF for 1/fcpa = 2. (thick line), for which 
A/|/ij = 0.474. The threshold for pair-breaking is given by the thiner line. The dashed line indicates the result obtained from 
the Bogoliubov formula. 



VIII. APPENDICES 
A. Vicinity of /i = 

We give the results of the first order expansion in the vinicity of /i — 0. All the specific numerical integrals coming 
in this calculation are elliptic integrals [29] , which can all be expressed in this case in terms of Euler's gamma function 
r(l/4) ~ 3.62561 = f. One finds: 

' (^)" 2 = B + C>t (29) 



aV2^A A \ A J " A 

where A = 4 7 r 1 / 2 /f 2 ~ 0.5393, B = 1/{2A) ~ 0.9270, C = 3n A/8 ~ 0.6354 and D = l/(nA) ~ 0.5902. In particular 
one finds that, for fi = 0, A = A with A /E F = 47r 1 / 3 /r 4 / 3 ~ 1.0518, with a corresponding value: 

/ 1 \ 8tt 2 / 3 1 /A \ 2 , , 



In the vicinity of a , we have from Eq.(29): 



7T/4 2 1 / 1 \ Rt/3 ,, 

A-A = — -0.4569m -, h— ) = -(1 + 7r 2 A 4 /4) — — ~ -0.6956 (31) 



B. COLLECTIVE MODE AND RESPONSE FUNCTION 



In order to obtain the response function various approaches can be used, which are actually all basically equivalent 
presentations of the standard self-consistent BCS theory. We will consider two of them. One is the Green's function, 
or diagrammatic, approach to the BCS theory. The other one is kinetic theory In the following subsection we will 
use the Green's function approach to calculate the response, and more precisely the condition under which the vertex 
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diverges. Then in the next subsection we will use the kinetic equation method to obtain the condition of divergence of 
the density response. Naturally each one of this method can be applied to the calculation of the full response, which 
we will give below. In both cases we will be rather sketchy since these methods are well-known in the literature, and 
we give only the following indications for completeness. 

1 . DIA GRA MM A TIC A PPR OA CH 

We will obtain the vertex T in the superfluid state by writing a Bethe-Salpeter equation [30]. This vertex is 
a generalized T- matrix. In the present case it depends only on the total momentum q and frequency w, as can 
be seen in the explicit setup of the Bethe-Salpeter equation, done below. More specifically we look for a pole of 
this vertex. We introduce the normal and anomalous Green's functions, G and F respectively, in the superfluid 
phase. It is more convenient to use the Euclidean form (after Wick transformation u) — ► ioj) for which we have 
G(iu,k) = -(iuj + £ k )/(u; 2 + E 2 k ) and F(iu,k) = -A/(uj 2 + E 2 k ). 

In order to obtain a Bethe-Salpeter equation, we have to introduce two vertices Tn and Ti2 (see Fig. 9). The first 
one Fn corresponds to the scattering of two atoms with opposite spins (with two in-going lines and two out-going 
ones), while the second one Ti2 has four out-going lines. This last one is non zero only in the superfluid state, while 
Tn exists also in the normal state. 




FIG. 9. Two vertices Tn(q) and Ti2(q) in superconductive state, p = (p, ui) and q — (Sl,q) are 4-momenta 

We obtain the Bethe-Salpeter equation by considering first the simplest process corresponding to an elementary 
interaction, depicted as the first term in the right-hand side of Fig. 10. Then we consider more complicated processes 
obtained by repeating this process, producing in this way ladder diagrams. In the intermediate steps we can introduce 
normal Green's functions, as in the second term in the right-hand side of Fig. 10 (they are depicted by lines with 
single arrows). But we can also introduce anomalous Green's functions, as in the third term in the right-hand side 
of Fig. 10 (they are depicted by lines with double arrows having opposite directions). On the other hand we can not 
have an intermediate step made of one normal propagator and an anomalous one because this would imply interaction 
between atoms having the same spin, which is forbidden for s-wave scattering by Pauli principle. Summing up all the 
processes after the first one lead to the equation represented diagrammatically in Fig. 10. In algebraic form it reads 

T 11 (q) = V-VT 11 (q) J^ G (p" + q/2)G(-p" + q/2)-VT 12 (q) J ^ F{p" + q/2)F{-p" + q/2) (32) 

where q is now a four- vector representing (w,q). If we want to take " renormalization" into account, as we have 
done in section II when going from Eq.(l) to Eq.(2), it is easy to show that we should replace GG in Eq.(32) by 
GG — Go Go, where Go(w,k) = [ioj — £fc] -1 is the vacuum Green's function. Simultaneously we have to replace the 
bare interaction V by the coupling constant g = Aira/m. From Eq.(32) it is natural to introduce now the elementary 
response functions Xij corresponding to the various bubbles appearing in Fig. 10. Specifically we define: 

-Xii(q) = \ + J-£fr[G(p + q/2)G(-p + q/2) - G (p)G (-p)} 

(33) 

Xu(q) = J A F(p + q/2)F(-p + q/2) 
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FIG. 10. The Bethe-Salpeter equation for the vertex Fu 



Then equation (32) takes the form: 



ri 2 (<?)xi2(<?) = 1 + rii(g)xii(g) 



(34) 



Now we have to derive a second Bethe-Salpeter equation for T\2(q). Graphically it has the form shown on Fig. 
11. A small difference is that the anomalous Green's functions F + (p) appearing in this equation have their arrow 
pointing outside, instead of inside as in Fig. 10. However for real A, they are simply related to the preceding ones by 
F+(p) = F s (-p). 

In algebraic form it reads: 



ria(g) = - 



<?r 12 (g) J |^ G[p" - q/2)G(-p" q/2) - gT^q) J |^ F+(p" + q/2)F+(-p" + q/2) 



(35) 



Renormalization requires again substitution of GG by GG — GqGo and the replacement of V by g = Aira/m. We see 
naturally appearing in this equation the above response functions Eq.(34), but for the 4- momentum —q, leading to 
the following form for this second equation: 



ri2(£?)xn(-£?) = r u (<7)xi2(-<?) 



(36) 



From Eq.(34) and Eq.(36) the vertices T\\ and Ti2 are immediately obtained. However we are only interested in the 
pole of these vertices, for which they diverge. This allows us to solve only the homogeneous equations corresponding 
to Eq.(34) and Eq.(36), leading to the following equation for the mode: 



Xii(<7)Xn(-9) = Xia(?) 



(37) 



where we have used the fact that, taking the explicit form for F{k) in Eq.(33), one has Xli{~ q) = Xi2{l)- This result 
has already been derived in Ref. [12]. 

In order to obtain a more convenient equation for the mode, we can perform the integration over the frequency 
variable in Eq.(33). This is easily done by closing the contour in the upper half- plane where the quantity to be 
integrated has two poles, located at iE + and iE—. Here we have introduced the convenient notation E± = Eh± q /2, 
and we will use similarly £± — £,k±q/2- The results are: 



■Xn(w,q) = J / 



d 3 k 

(2tt) 3 



(E+ + E_)(E+E_ + e+C-) + iuKE+Z- + g-Cf) 1 
E+E_ [(E+ + E_) 2 + iv 2 ] E k 



(38) 



, 1 f d 3 k A 2 E + + E- 

Xl2( "' q) = 2 J (2.y E+ E- [{E ++ E-f + ^] m 

where, in Eq.(38), we have made used of the gap equation Eq.(2) to get rid of the term 1/g in Eq.(33). Finally we 
have to go back to ordinary frequencies by the " inverse" Wick transformation iuj — > uj to recover the results of section 
IV for the pole corresponding to the collective mode. 
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FIG. 11. The Bethe-Salppeter equation for Ti2 



2. KINETIC EQUATION 

This kinetic equation approach is the generalization of the well known Landau's kinetic equation for Fermi liquids 
to the superfluid phase described by BCS theory. It provides a dynamic self-consistent mean-field theory which is 
gauge invariant and satisfies from the start conservation laws. In particular, in the appropriate limit, it reduces to 
hydrodynamics as it should. It has been derived from the BCS formalism by Stephen [31], generalized by Betbeder- 
Matibet and Nozieres [32] to include Fermi liquid effects and later used by Woffle [33] in the case of superfluid 3 He. Its 
application to cold gases has been pointed out by Randeria [16]. A specific application close to the present derivation 
where details can be found is Ref. [34]. Here we use only the collisionless version of this equation. 

To obtain easily this equation in a compact form it is convenient to write the standard BCS Hamiltonian, which 
rules the dynamics of the system: 



hbcs = 4 Q c fcQ + X)( A 4t c -fci + h - c -) 



k <v 



with Nambu spinor notations as (within an additive constant): 

Hbcs = ^2 c h<- 

kXfi 



(40) 



(41) 



where A, fi = 1 or 2, and by definition Cf. i = t and c^- 



This corresponds to a unitary transformation and 



implies the anticommutation relations [ck\,c\, ] = 6k,k' and [ck\,Ck> M ] = 0. The 2x2 matrix ekx^ is given by: 



(-k 



(42) 



and tr indicates the transposition. 

We consider now a small perturbation 5U with frequency uj and wavevector q, coupled only to density fluctuations 
for simplicity. This gives rise to a corresponding fluctuating part SA of the order parameter, with the same frequency 
and wavevector. These perturbations give an additional contribution to the Hamiltonian: 



5H= ^ 



/ 4 C k+q/2,\ 06 V C k-q/2.,fi 
kXfi 



(43) 



with Sin — —5^22 — SU and <5e2i = (#£12)* = <5A. The time evolution of the single particle density matrix: 

(?WW = (4A c fe>) ( 44 ) 

can then be obtained from Heisenberg equations of motion and commutation relations. One finds that the fluctuating 
part of this density matrix 8hk — nk-q/2, k+q/2 satisfies the following kinetic equation: 



17 



ujSh k = Sn k e° k+q/2 - e° k _ q/2 Sn k + n° k _ q/2 Se k - 5e k n° k+q/2 



(45) 



Here e k is the equilibrium value of the matrix i k given by Eq.(42), A taking its equilibrium value, and h k is related 
to the corresponding T = equilibrium value n kk , of n kk > given by h kk , — n k S kk >, and it is given by: 



E k 



(46) 



Naturally the order parameter is self-consistently related to the off-diagonal part of h k k > , both for the equilibrium as 
well as for the fluctuating part. These equations are conveniently combined to eliminate the interaction strength V 
and yield the self-consistent relations: 



SA 
2E~ k 







SA* 
2El 







(47) 



which allows to close the set of equations. In this way the response to the perturbation SU can be calculated, and its 
pole provides the mode dispersion relation. In the present case, since we are only interested in this pole and not in 
the response function, we may set from the start SU = and the equation for the mode dispersion relation will just 
be obtained from the compatibility condition for Eq.(47.) The algebra is simplified if one prcmultiplies the matrix 
equation (45) by U k _ q / 2 and postmultiplies it by U k+q / 2 where: 



U k 



Uk v k 
-v k u k 



(48) 



where u k = [(1 + £ k /E k )/2] 1 / 2 and v k = [(1 — £ k /E k )/2} 1 / 2 are the well-known coefficients of the Bogoliubov-Valatin 
transformation, which diagonalizes both e k and n k . This leaves only diagonal elements for the matrix equation which 
is then easily solved. We skip this algebra and give only the result for the off-diagonal elements of Sn k (oj, q): 



(Sn k ) 2 i(uJ,q) = {Sh k )l 2 (-<jj,q) 



2 uj 2 - {E++E-Y 



SA\^{^ + ^) + {E + +E_){l + ^\ 

A2 

- SA*(E+ + E_) 



(49) 



E-E+ 



Intermediate steps can be found in Ref. [34] . It is easy to check that this equation gives the proper result in the static 
and uniform limit ui = 0, q = 0. Finally it is worth noticing that Eq.(45) yields automatically the mass conservation 
equation, by taking the trace and summing over k. 



3. RESPONSE FUNCTION 



Here we give for completeness the explicit expression of the density-density response function x( ( 7j u ') = Sn/SU, 
where Sn = J2 k (^ n k) n is the single spin state density fluctuation produced by the perturbation SU. In addition to 
the integrals introduced in section IV we set: 



J(w, q) = [ k 2 dk [ du—t 
Jo Jo E + 



(Z++Z-)(E + +E_) 



E_[(E+ + E-) 2 - lo 2 } 



(50) 



hi(u),q) - J 2 2(^,g) 
2A 2 



(51) 



I"(u,q) - rk 2 dk ( du^±- 
Jo Jo E + 



E- (E+ + E-) 2 



Then the response function is given by: 
1 



xO?,w) = - 



4tt 2 



I"(u J ,q)-A 2 



\ (E+ + EJ) E+E- - + A 2 

-uj 2 

hi(u;,q)I 2 (Lj,q) +oj 2 I 22 (u,q)r 2 (LJ,q) - 2uj 2 Ii 2 (uj, q)I(w, q)I'{u, q) 



hi (oj, q)I 22 {uj, q) - uJ 2 lf 2 (oj, q) 



(52) 



(53) 
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The first term, proportional to I" is the simple and naive BCS result one would obtain by applying the standard linear 
response formula to the BCS ground state. Its poles exhibit, at low momenta, the typical gap 2A of BCS theory. 
The second term is the contribution coming from the dynamics of the order parameter, and plays a crucial role in 
providing the "phonon" dynamics absent in BCS theory. In particular this term restores the /-sum rule, which results 
from particle conservation and which is violated by simple BCS theory, and its poles yield naturally the dispersion 
relation Eq.(14) of the collective mode. Note that the imaginary part of the response function Eq.(53) is proportional 
to the dynamic structure factor whose behaviour along the BEC-BCS crossover has been recently investigated in [22]. 



C. VICINITY OF THE CONTINUUM ON THE BCS WEAK COUPLING SIDE 



One can make, in the limit ui — ► 2A where the collective mode is very near the continuum threshold, an analytical 
study of Eq.(14). In this limit the denominator (E + + EJ) 2 — lo 2 in the integrals In,I\2 and 7 2 2 is nearly zero 
when £ + and £_ goes to zero. This condition on £ + and is equivalent to require that u is small and k is near 
kth — \/2rn\i — <7 2 /4. However a closer examination shows that the numerators in I\2 and J22 go also to zero in 
this case, and as a result these two integrals have actually a finite limit when lo — > 2A. On the other hand there is 
indeed a divergence which appears in In- Since, from Eq.(14), In has to be finite when the dispersion relation is 
satisfied, this divergence has to be cancelled in some way. Taking for example Eq.(18) for In, one sees that only the 
lo 2 term has a divergence. It is easily seen that, within a constant prefactor, the divergent term of In behaves as 
(1/q) ln[l/(2A — ui)]. In order to compensate for the logarithmic divergence, q has to grow like ln[l/(2A — lo)] when 
co — > 2A. This means that the dispersion relation Lo(q) approaches the limit lo = 2A with an exponential behaviour. 
This result is completely consistent with numerical results. 

Naturally this result holds only as long as kth exists. This implies as a necessary condition that q < 2\f2mjl. One 
can consider directly the case q = 2yJ2mfi. It is easily seen that the qualitative situation is not changed, i.e. In is 
divergent for lo — 2A whereas the other integrals are finite. This implies again that one must have lo < 2A in order 
to satisfy the mode dispersion relation. On the other hand the quantitative situation is slightly different because the 
divergence is produced by contributions to the integral coming from the vicinity of k = 0. Then, when one goes in 
the domain q > 2 v / 2m/z, the situation changes qualitatively because as we have seen in section VA the minimum of 
E + + E_ does not correspond anymore to the condition £ + = = 0, but rather to £,+ =£,-= q 2 /(8m) — /1, obtained 
for k — 0. Because of the phase space factor k 2 in all the integrals, the three integrals / are finite at the continuum 
threshold. As it is seen in section VC, the collective mode does indeed merge in the continuum in the weak coupling 
regime, but this occurs always for q > 2y/2mJI. However when one goes toward weaker coupling, this merging point 
is increasingly close to the q = 2y/2mjl point. 



D. FIRST ORDER EXPANSION ON THE BEC SIDE 



We want to expand the dispersion relation Eq.(14) on the BEC side to first order in (A/\fi\) 2 . This means going 
one step further than what we have done in the preceding calculation in section VD. This turns out to be not so 
complicated analytically. Indeed the only different term between In and / 22 , namely the A 2 term in Eq.(16) and 
(17), gives a negligible A 4 contribution in the product /11I22 coming in the dispersion relation Eq.(14). This implies 
that we can omit this A 2 term in Eq.(16) and (17) and take In = h2- Hence the dispersion relation simplifies again 
into In + w/i2 = as precedingly. Moreover we have to first order in A 2 , E± = £± + A 2 /2£±. This implies that we 



have at this order, for the numerators coming in In and J12, E+£- + £,+E_ 
relation becomes: 



f k 2 dk f 
Jo Uo 



du(l 



E-\-E— E-\- -\- E- 



2lo 



1 

E 



E + E- + so that the dispersion 



= (54) 



Here we use reduced units where all the energies are expressed in units of and all wavevectors in units of k n = 
\j2m\n\, with the notations £1 = lo/\/i\ and Q — q/ko. All the quantities in reduced units are indicated by a bar over 
the corresponding symbol. In terms of these reduced units, a first order expansion leads to: 



A 2 [°°- k 2 

/ dk =7- 

4 ./n 1 + k 2 



du 



10 ' Jo 

After performing the integrations one obtains: 



( £l l + fc 2 ^+ (1 + k 2 ) 2 



(55) 



19 



Q 2 A 2 16- Q 2 

"-"TieTQ* (56) 



equivalent to: 

O 2 = n 2 

4 16 + Q 2 V 2 

up to first order in A 2 . 



«= Q ^^ + m 2 tan 
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